A New SUSY mass reconstruction method at the CERN LHC 



K. Kawagoe, 1 M.M. Nojiri, 2 and G. Polesello 3 

1 Department of Physics, Kobe University, Kobe 657-8501, Japan 
2 YITP, Kyoto University, Kyoto 606-8502, Japan 
3 CERN, CH-1211 Geneva 23, Switzerland and INFN, 
Sezione di Pavia, Via Bassi 6, 21100 Pavia, Italy 
(Dated: February 2, 2008) 

We propose a new mass reconstruction technique for SUSY processes at the LHC. The idea 
is to completely solve the kinematics of the SUSY cascade decay by using the assumption that 
the selected events satisfy the same mass shell conditions of the sparticles involved in the cascade 
decay. Using this technique, we study the measurement of the mass of the bottom squarks in the 
cascade decay of the gluino. Based on the final state including two high pr leptons and two b-jets, 
we investigate different possible approaches to the mass reconstruction of the gluino and the two 
bottom squarks. In particular we evaluate the performance of different algorithms in discriminating 
two bottom squark states with a mass difference as low as 5%. 

PACS numbers: 



I. INTRODUCTION 

Supersymmetry is one of the most attractive models beyond the Standard Model (SM) of elementary particle 
physics. The superpartners of SM particles (sparticles) might have masses of the order of the TeV, and experiments 
at the Large Hadron Collider (LHC) should be able to detect such particles up to masses of 2-3 TeV 0,0- 

The pattern of the sparticle mass spectrum depends on the SUSY-breaking mechanism, which might depend on 
gravity, space-time structure, or unknown interactions. The unraveling of such mechanism through the determination 
of the sparticles' masses is therefore one of the most important physics targets of future collider experiments. 

The potential of the LHC for SUSY mass determination has been studied in detail in the past decade. The 
most promising method involves the study of the endpoints in the distributions of invariant masses among the visible 
sparticle decay products. Information on the masses involved in the cascade decay can be extracted from the endpoint 
values if the decay distributions are dominated by a single cascade decay chain. Studies based on the endpoints are 
documented in 0, 11 |jj |1 . 

In this paper we explore a new method for reconstructing SUSY events. This method does not rely only on events 
near the endpoint. Instead, one kinematically solves the neutralino momenta and masses of heavier sparticles using 
measured jet and lepton momenta and optionally a few mass inputs. 

We concentrate in this paper on the measurement of the mass of the bottom squarks (sbottoms) through the cascade 
decay: 

~g^bb 2 -> x° 2 hb 2 - lb x b 2 t 2 -> (1) 

We address this decay, rather than the equivalent cascade involving a generic qr, for various orders of reasons. First 
of all, the physics of the third generation squarks and leptons is particularly important in disentangling the pattern 
of SUSY breaking. In the MSUGRA model, hi, is lighter than b r due to the RGE running by top Yukawa coupling. 
The mass of the lighter sbottom state b\ is further reduced by the mixing of the left and right b states. The sbottom 
sector is thus sensitive to fundamental parameters of the theory such as the trilinear couplings and tan (3 otherwise 
of difficult access at the LHC. The third generation sparticle masses are also important parameters for B and Higgs 
physics. Second, the chain in Eq. Q involves two b quarks, which can be tagged in the detector. The problem of 
correctly identifying the jets contributing to the interesting decay is thus made much simpler. Finally, as we will 
discuss in detail below, both sbottom states b\ and b 2 yield the decay chain of Eq. (JTJ. The performance and the 
robustness of the reconstruction algorithms can be benchmarked against the ability in disentangling the two states. 

Five sparticles are involved in the cascade decay Eq. (JJ, therefore one can write five mass shell conditions among 
the leptons and quarks in the final decay products. 
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m S = ( p x? +Wi +PI2 +Pbi +Pb 2 ) 2 - (2) 

For a bb££ event, the equations contain the 4 unknown degrees of freedom of the x? momentum. Each event therefore 
describes a 4-dimensional hyper-surface in a 5-dimensional mass parameter space, and the hyper-surface differs event 
by event. From the purely mathematical point of view 5 events would be enough to determine a discrete set of 
solutions for the masses of the involved sparticles, and the probabilistic discussion of the following could be easily 
developed in a 5-dimensional space. In order to illustrate the method in a more transparent way, we will develop the 
argument by assuming that the masses of x,2> £, and Xi are known. This is a reasonable assumption at the LHC, 
where it has been shown that a detailed study of the lepton-lepton system from the X2 decay can be used to precisely 
constrain these masses |8j- In this case, each event corresponds to a different line in the (nig, mg) plane, therefore two 
events are enough to solve the gluino and sbottom masses altogether. 

We call this technique the "mass relation method" , because here one uses the fact that sparticle masses are common 
for all events which go though the same cascade decay chain. Note that the events need not to be close to the endpoint 
of the decay distribution, but they are still relevant to the mass determination. This means that one can use the mass 
relation method even if the number of signal events is small. 

The purpose of this paper is to explore in detail the implications of the mass relation method which is only sketched 
for signal events in 0. In that report a measurement of the gluino and sbottom masses was obtained from the peak 
of the distribution of the solutions for all possible event pairs, assuming that the X2> £ and Xi masses are known. In 
this paper we extend the previous analysis to take SUSY backgrounds into account in the distributions. 

We note that b in Eq. Q could be either bi or b 2 . The decay was studied in Q by using the endpoint method, where 
the possibility of distinguishing the two sbottom states b\ and b 2 was studied for a case where the mass difference 
was approximately 5% of the sbottom mass. The conclusion showed that even with very large integrated luminosity 
the result is at best marginal, and crucially depends on the ability of the experimenters to model to a very high level 
of detail the response function of the detector to 6-jets. This is probably an inescapable conclusion, given that the 
resolution in calorimetric measurement is comparable to the mass splitting one wants to measure. It is however worth 
studying alternative reconstruction methods which make a better use of all the information available for the selected 
events. 

In a subsequent paper |lOj ]. a similar process is analyzed by constructing an approximate event-likelihood function 
for signal events, based on a Bayesian statistical approach. In this paper we construct a different approximation for 
the likelihood function. We apply this function to the analysis of a simulated data set, and we study how well b\ and 
62 can be reconstructed for a sample model point and its variants. 

This paper is organized as follows. In section [HI we introduce the detailed SUSY models addressed in our study, 
and briefly present the simulation procedure. We then solve in section ITTT1 analytically the kinematics for the process 
Eq. , arriving to a compact expression for the masses of g and b as the function of b and £ momentum and discuss 
the typical solutions. In section [V] we discuss the event pair analysis on our sample SUSY model, where rrig and 
are computed from all event pairings in the selected samples. Gluino and sbottom masses consistent with the 
input values are reconstructed, but we are also forced to artificially select one of the multiple solutions obtained 
from solving for the masses the coupled quadratic equations for each event pair. In section IVTl we therefore define 
an approximate likelihood function built using all the events in the sample, and we describe an analysis based on 
this function. The masses reconstructed from the likelihood analysis are in agreement with the input values and 
the method automatically takes care of the issue of multiple solutions without artificial selection. We also study 
the possibility to reconstruct the b 2 in our sample model. Section IVTT1 is devoted to discussions and comments. We 
especially discuss the theoretical relevance of being able to perform detailed measurements in the third generation 
squark sector. 



II. MODEL POINTS AND SIMULATIONS 

We choose for this study the model point SPSla as defined in Ref. and its variants. SPSIa has a significant 
production cross-section for the chain of Eq. QJ, and has already been the subject of detailed LHC studies 0. 
The model is defined in the mSUGRA scenario by the parameters tan/5 — 10, m = 100 GeV, M = 250 GeV, 
Aq = —100 GeV, and [i > 0. In order to evaluate the performance of the analysis for different values of the splitting 
of the two sbottom states, we also study the points where tan/3 = 15 and 20, keeping other GUT scale parameters 
same as those for SPSIa. For the additional points, the g, X2, Xi and £ masses are almost the same as for SPSIa, 
while mg is reduced because the left right mixing of b increases proportional to \x tan (3. 

The masses and decay branching ratios are calculated with the ISASUSY code The mg changes by 3% from 

tan/3 = 10 to tan/3 = 20. We will see later that such a difference in the 61 masses should be measurable at the LHC 
if systematic effects are kept under control. 
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Another drastic effect in Table ITT| is the decrease of Br(x\ — ¥ f£) = e, fi). This is because the decay width for 
X2 — * tt increases as the f mixing angle increases. For tan (3 = 20, the branching ratio into bb££ mode is reduced by 
a factor of 5 yielding a rather small number of accepted events. 

We also list in Table ITT1 the decay branching ratios relevant to the decay chain of our interest. The ratio of the 
decay branching ratios — Br(g ~ * 62 — ► i£i)/Br(g — > b\ — > x§) varies from 0.26 to 0.14 when tan/3 changes from 
10 to 20. This tan/3 dependence comes from the reduction of b\ mass and the increase of the left-right mixing of b, 
9, defined as l x = cos 9b L + sm0b R , from 9 = 0.49 (tan/3 = 10) to 0.61 (tan/3 = 20). Namely, 

• T(g — > 61) increases as is reduced. 

• r(& 2 — > X2) increases as 9 increases. However, T(b 2 — > Xi") and ^(62 - * ?i) also increase, therefore Br(b 2 — > X2) 
stays more or less same. 

• is reduced and T(6i — > t\W) is kinematically suppressed as tan/3 increases. 

The branching ratio is in principle important information to determine the mixing angle 9 but one needs to know the 
relevant sparticle mass spectrum to utilize it. 

Note that the ratio is small and the mass splitting between b\ and b 2 is within 10% over the parameters we 
study. The b\ is therefore significant background for the 62 search at the LHC. Nevertheless, a hint for the decay of 
g — > f>2 might still be observable in the bb££ distribution, given an excellent control of the detector response to 6-jets, 
as we will see in section IVT1 

The SUSY events are generated using theHERWIG 6.4 00 event generator. The produced events are passed 
through the ATLFAST jl5j detector simulator, which parameterizes the response of the ATLAS detector. In particular, 
we use the parametrization for 6-tagging efficiency corresponding to the expected high luminosity performance of the 
ATLAS detector. While the performance is a function of the pr of the jets, a typical performance figure is Cb = 0.5 
for the b tagging efficiency for a rejection of 100 on light quark jets. 

For each point we have generated 1.5 x 10 7 events which approximately correspond to an integrated luminosity of 
300 fb" 1 . 

The following cuts are applied in order to select signal events: 

• M e ff > 600 GeV and i?™ ISS > 0.2M c g, where £J™ ISS is the missing transverse energy and M c g is the scalar sum 
of the missing transverse energy and the transverse momenta of the four hardest jets, 

• at least 3 jets with p T1 > 150 GeV, p T 2 > 100 GeV and p T3 > 50 GeV, 

• exactly two ^-tagged jets with pr > 50 GeV, 

• exactly two opposite-sign isolated same-flavor (OSSF) leptons with pti x > 20 GeV, pre 2 > 10 GeV, with an 
invariant mass 40 GeV< m« < 78 GeV. 

The isolation criterion consists in requiring a transverse energy deposition in the calorimeters smaller than 10 GeV 
in a (77, 4>) cone of radius 0.2 around the lepton direction, where 77 is the pseudorapidity of the lepton and cf> the angle 
in the plane transverse to the beam. A detailed discussion of the Standard Model backgrounds after these cuts is given 
in ■ The authors show that the Standard Model background is negligible in comparison to the SUSY background, 
consisting of SUSY events not including the decay chain of Eq. . 

In order to perform the reconstruction, we need to identify the position of each of the two &-jets and of each of the 
two leptons in the decay chain. In the following analysis, we assume that the 6-jet with larger px originates from the b 
decay. The assignment is optimal for SPSla, because — m^o 3> wig — m^. At the same time one can fix the lepton 

assignment so that the higher (lower) pt lepton ^high (^low) comes from I to increase the possibility to pick up the 
correct lepton assignment. If we roughly know the masses of £,Xi an( i X2j ^ ^ s eas Y to determine which assignment 
is optimal For the likelihood analysis on SPSla, we use £hi g h as £\ in Eq. (JIJ. 

III. SOLUTIONS OF THE DECAY KINEMATICS 
A. Formula of the decay kinematics 

It is straightforward to solve the decay process in Eq. (|TJ . We first note that the b cascade decay can be written as 
a function of momenta of b\ , £\ and £2 , because the four mass shell conditions can be used to eliminate the unknown 
Xi four momentum. 
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m in 


x i 


10 


595.2 


491.9 (524.6) 


176.8 


143.0 


96.0 


15 


595.2 


485.3 (526.9) 


177.9 


143.0 


96.5 


20 


595.2 


478.7 (531.2) 


178.5 


143.1 


96.7 



TABLE I: Masses of the relevant sparticles for the three studied points, in GeV. 



tan P 


Br(g -> fei (2) ) 


Br(b 1(2) -> x§) 


Br(x° 2 - £) 


10 


8.9 (4.9) x2 


35.8 (16.7) 


3.16x4 


15 


9.8 (4.6) x2 


37.8 (15.9) 


1.26x4 


20 


10.8 (4.1) x2 


38.9 (13.2) 


0.61x4 



TABLE II: Branching fractions for the decays used in the analysis, in percent. 



To systematically solve the system, we expand the Xi momentum with the observed momenta of b±, t\, and l 2 : 



p^o = ap^ + bp£ 2 + cp bl . 



(3) 



The expansion is possible if p^ l; pi 2 and p^ are independent from one another. The three on-shell conditions may be 
then rewritten as 
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(5) 



The three parameters, a, b and c are solved as functions of E^o if det M. ^ 0. By using the on-shell condition of xi, 
we obtain the quadratic equation of E^o as 
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E x1 




where Aij — [x;] T Al x [xj] — 8ij(i — 2)/m? and definition of xi is given in the appendix 



T M-^-S ij {i-2)/ml Bi 
As the b decay kinematics is solved, we now use the on-shell condition of g — > bb 2 to obtain 

Qumi + 2Q 12 m\m\ + Q 22 m\ + 2Q x m\ + 2Q 2 m\ + Q = 0. 



(6) 



(7) 



Q's are functions of the momenta of the leptons and b quarks, m^o, m x°: an d m i- The expressions for Q's are shown 
in Appendix. 

Mathematically, when there are two independent events, we have two independent equations of the form Eq. J7J. 
The equation is quadratic for m~ and m? and can be analytically solved 1 . There are up to four solutions for an event 
pair, and one of them must coincide with the true solution. 



Note that when some of the momenta are parallel, JV[ cannot be inverted. We ignore the possibility, because experimentally we always 
require the isolation of leptons from jets. 
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IV. ENDPOINT ANALYSIS AND MASS RELATION METHOD 

The decay chain g — ► bb 2 — > X^i&S — > Ib^l-i — > Xi&iWi^2 has already been studied in detail in Refs. 0,0] > using 
the endpoint method shown in Ref. pj. In this method, the masses of the %5 an d X2 are assumed to be known from 
the study of the kinematic edges in the decay. For events near the edge of the mu distribution the Xi is essentially 
at rest, and the momentum of the X2 can be approximated with the relation: 



(8) 



This formula is correct only at the endpoint of the three body decay X2 ~ * Xi^> but is nevertheless approximately 
correct near the edge of ££ —* Uxi a * SPSla. The sbottom mass can then be calculated by building the invariant 
mass of the approximate X2 obtained by using Eq. (JHJ with the leading 6-jet in the event. The parton-level result is 
shown in Fig.^a), where we plot the difference between the reconstructed gluino mass and the reconstructed sbottom 
mass, to minimize the smearing introduced by the approximation of Eq. |JSJ. The &i peak is reconstructed correctly. 
The bump on the left originates from the events from g — > 62, a nd is centered on the position nig 



mr 



=70.6 GeV. 

Even at parton level, with no experimental smearing the bump is not well separated from the b\ peak at 103 GeV. 
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a) b) 

FIG. 1: The distribution of m g — calculated using a) parton level b momentum and using the approximate relation Eq. JSJ 
and b) by solving Eq. For a) mu > 65 GeV is required. For b), the two solutions for the input gluino mass of 

nig = 595 GeV are plotted in a two-dimensional plane, for both correct and wrong lepton assignment. 

Unlike the relation Eq. (JSJ), the formula Eq. Q is exact. For each event we have in this case two possible lepton 
assignments, and for each lepton assignment, given an input value for the gluino mass, two solutions for the sbottom 
mass from the quadratic equation Eq. Q. We show in Fig. ^ the smaller solution versus the larger one of Eq. (JJJ for 
both of the lepton assignments. For the gluino mass the nominal value is assumed. Unlike in Fig.^a) two peaks of 
b\ and 62 a re clearly separated in the plot. Furthermore the number of the available events for the mass fit are now 
increased by factor of 2 compared to the endpoint analysis because there are no constraint on the value of mu when 
using Eq. (J7J. 

The advantage of switching to the exact solution for the event kinematics is clearly demonstrated by the plots. In 
order to evaluate if the heavier sbottom state will be detectable we need to perform the study taking into account the 
smearing induced by the fact that &-jets rather than partons are measured in the detector. 



V. EVENT PAIR ANALYSIS 



As discussed in section ITTT1 each event can be represented as a quadratic equation in the m|, m~ variables. By 
taking two events, we have a system of two equations in two unknowns which can be solved analytically. This yields 
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tan/3 = 10 


tan P = 15 


tan/3 = 20 


m- a - m- b 


591.9 
98.9 


593.1 
105.1 


585.1 
111.6 



TABLE III: Fit results of the gluino and and sbottom masses in the event pair analysis. 

up to four values for the squark and gluino masses. In the following, We start from the selected bbU events, and 
we build all the possible event pairs. In order to minimize the combinatorial backgrounds we use the pairings which 
satisfy the following conditions 2 : 

• Eq. 10) has solution for only one of the two possible lepton assignments. 

• For the selected lepton assignment the resulting quartic equation in m| has only two solutions, and the difference 
of the gluino masses for the two solutions is more than 100 GeV. The smaller gluino mass solution is chosen. 

The rrig distributions for the OSSFxOSSF events pairs are shown in the histograms on the upper line of Fig. O A 
significant SUSY background, also shown in Fig. [3 is still present in the sample. This background can be estimated 
from the data themselves by using the bbil events with an opposite sign opposite flavor (OSOF) lepton pair (i.e. 
U = e fj.^). To this purpose we produce mass distributions for the three types of event pairs: 

1. two OSSF lepton events (OSSFxOSSF), 

2. an OSSF lepton event and an OSOF lepton event (OSSFxOSOF), 

3. two OSOF lepton events (OSOF x OSOF). 

The background-subtracted distribution can be then obtained as the combination of the three distributions: 
OSSFxOSSF-OSSFxOSOF+OSOFxOSOF. 

The total number of event pairs is of course much larger than the number of events. Due to the selection criteria 
imposed to minimize the combinatorial backgrounds, some of the events are not used at all to make event pairs, while 
events which are used at least once are used O(10) times on average. The three histograms show peaks corresponding 
to the input value for the gluino mass even before the background subtraction. The green and blue histograms show 
the OSSFxOSOF and OSOFxOSOF distributions, respectively. The distributions after the background subtraction 
are shown in the histograms on the lower line of Fig. |2 The peak position and its error obtained by a Gaussian fit 
to the distribution are listed in Table ITTT1 We also show in the histogram the gluino mass distribution for event pairs 
where at least one of the events comes from 62 decay. The &2 contribution is significantly smaller compared with the 
distribution of b± pair and does not affect the gluino mass fit. One can also look into the distribution of nig — and 
estimating the value of this observable by performing a Gaussian fit on the observed peak. The result of the fit is also 
shown in Table UTTl 

The statistical error of the gluino mass measurement can be evaluated by performing the analysis on a set of 
statistically independent experiments. To perform the evaluation within a reasonable CPU budget, we generated 
a set of events where the b\ is forced to decay with 100% BR into the desired decay chain, for the SPSla Point 
(tan/3 = 10). The generated statistics corresponds to 30 experiments with an integrated luminosity of 300 fb" 1 each. 
By performing the analysis on the 30 experiments, we find that the spread for the measured gluino mass is 1.6 GeV 
for the benchmark integrated luminosity of 300 fb _1 . By construction, the presence of the combinatorial background 
is not considered in the error analysis. No large degradation of the resolution is expected once this effect is correctly 
taken into account. 

Once the gluino mass is fixed by the analysis shown above, Eq. Q can be solved for each event for the sbottom 
mass, giving as input the central measured value for the gluino mass. Following this procedure, in Fig. [21 we plot the 
distribution of the smaller sbottom mass solution mz(vam) for both OSSF (signal) and OSOF (background) lepton 
pair events (top histograms). We show in the bottom line the mass distributions after the background subtraction, 
where the red regions show the contribution of bi. The peak positions are evaluated by a Gaussian fit, and listed in 
Table ITVl Note that the total number of the signal b\ event is smaller by factor of 4 for tan/3 = 20 compared with 
tan/3 = 10, but the mass peak is seen very clearly. The m^min) peak and m^(input) are in good agreement, and it 
is not so for the larger solution m^max). 



2 Note the selections are rather phenomenological and they may introduce some bias to the reconstructed sparticle masses. We find 
however the obtained peak positions are consistent with the input masses in this study. 
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FIG. 2: The gluino mass distributions for tan/3 = 10 (left), tan/3 = 15 (center) and tan/3 = 20 (right) with the event pair 
analysis. The open, green, and blue histograms in the top figures are for OSSFxOSSF, OSSFxOSOF, and OSOFxOSOF event 
pairs, respectively. The open histograms in the bottom figures shows the mass distributions after background subtraction. The 
contributions of &2 are shown by red histograms. 





tan/3 = 10 


tan/3 = 15 


tan /3 = 20 


mj(true) 


491 


485.3 


478.8 


mi (min) 


492. 1± 1.2 


487.7± 2.2 


474.3± 2.4 


mj(rnax) 


504.5± 1.0 


502.9± 1.7 


495. 1± 2.4 



TABLE IV: Fit results of the sbottom mass with a fixed gluino mass (m§ = 595 GeV). 



The peak positions of the distribution of the events originated from 62 decay are also consistent with 62 masses. It 
is however evident from the plots that in the real experiment it will not be possible to claim the presence of the second 
peak. The existence of 62 can be established only after understanding &-jet smearing and b\ distribution correctly. 



VI. LIKELIHOOD ANALYSIS 



The analysis shown in the previous section only uses the events in pairs to evaluate the values of the squark and 
gluino masses. A more efficient use of the available event statistics can be achieved by using all of the events at the 
same time and finding the [rrig, mi) pair for which the combined probability of all events is highest. We define in this 
section an approximate likelihood function for [mg,mi), and we then apply it to detection of 62- 



A. Construction of the likelihood function 



From Eq. JJJ, each event is represented as a curve in the (m g , m^) plane. The coefficients of the curve are a function 
of the four momenta of the detected partons. The partons are measured as jets in the detector, which smears the 
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FIG. 3: The distributions for tan/3 = 10 (left), 15 (center) and 20 (right) with a fixed gluino mass (mj = 595 GeV). The 
open and green histograms in the top figures show the distributions of OSSF and OSOF lepton events, respectively. The mass 
distributions after background subtraction are shown in the bottom figures, where the contribution of the 62 events are shown 
by the red regions. 



parton according to a smearing function. It depends on the detector performance and on the algorithms used to 
cluster the energy deposition in the detector. 

Due to these experimental effects, in a frequentist approach, from the measured quadratic form for an event we 
can build a confidence belt in the (m^nij) plane. This should be built using a Neyman construction |17| . from 
the probability distribution for the measured values of (m^m^) as a function of the input values for the same two 
variables. In order to build this function, the crucial ingredient is the distribution of the measured 6-jet momenta as 
a function of the 6-parton momenta. Evaluating this distribution is outside the scope of this work, as it requires a 
detailed simulation of the detector response, which will need to be validated on real data using calibration samples of 
6-jets in the detector. Even assuming an approximate form for the 6-jet response, the proper Neyman construction 
for each event is a very computing-intensive calculation. 

We recast equation Eq. Q in the form 

f{mg,m i ,p 1 ,p 2 ) = 

where p\ and p 2 are the momenta of the two 6-jets. We build an approximate probability density function according 
to the formula: 

£(m s ,m~ b ) = J dp[ J dp' 2 (.(p Y : p'i)e{p 2 : p' 2 )S(f(m g ,m~ b ,p / 1 ,p 2 )) (9) 

where e(p : p') is the probability to measure a momentum p' for a 6-jet, given a 6-parton with momentum p. 

In the equation we did not include the possibility of lepton momentum mis-measurement, which has an almost 
negligible effect, as compared to the smearing of 6-partons. As a further simplification, we assume that the jet 
direction is not modified by the measurement and we use for e(p : p') a gaussian distribution, with a width a 
corresponding to the parameterized jet smearing used in the fast simulation program. 



cr/E = 0.5/V£(GeV) +0.03 < 3.0) 
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a/E = 1.0/y/E(GeV) +0.07 > 3.0) (10) 

The gaussian smearing is not very good approximation for for 6-jets for which in many cases the semi-leptonic decay of 
the 6-quark results in jets containing an unmeasured neutrino. The approximate function takes however into account 
the dominant part of the jet smearing and can be used to demonstrate the method. Note also that in Eq. [5J we 
define our C using e(pi,p[) where p\ is measured 6-jet momentum. The function C would correspond to the actual 
probability function only if the jet response were gaussian. A detailed experimental simulation will be needed in order 
to assess the validity of the obtained results in the real experimental situation. 

We now show log£ in the (m g — mg, nig) plane for a few events where the bb££ events originates from the cascade 
decay of Eq. (JTJ at SPSla. We calculate C using the following procedure. For each event in our sample, characterized 
by a (pi,P2) pah of measured momenta for the 6-jets, we generate Monte Carlo events where the two 6-jets have 
momenta (p\, p' 2 ), where p\ and p' 2 are randomly generated according to the function e{pi : p'{) x e(p 2 : p' 2 ). Each 
generated event corresponds to a curve in the (m g , mt) plane which satisfies the equation / (nig,m b ,pi, P2) = 0. We 
histogram of the number of curves that go through each bin of a 1 x 1 GeV grid in the (mg,mg) plane, for n Monte 
Carlo events, normalized by dividing the bin contents by n. In the limit n — > 00, this corresponds to 

£(m s ,m s + Amg;mg,mg + Amg) = J Dp[ J Dp 2 e(p 1 : p[)e(p 2 : p' 2 )9(p' l7 p 27 m g7 m g + Anig;m i ,m b + Amg) (11) 

where 0(p' 1 ,p' 2 , nig, m g + Amg; mg, mg + Anij) is 1 when the solution of Eq. Q for the two &-jet momenta p[ and p' 2 
goes through (nig, nig + Am^; mg, mg + Amg) and otherwise 0. We take n = 10000 for our calculations. 
In Fig. 0] we plot 

Alog£ = \og(C(m g ,m g + Am g , m~ b ,m~ b + Am~ b ) + c) — log(£(min)). (12) 

where c = 0.001 is a constant cutoff factor, which is needed as for each event we generate only a finite number of Monte 
Carlo experiments, and therefore some bins can have zero hits. The shape of the probability density distribution is 
different event by event, as it depends on the event kinematics. For a few events the density distribution is parallel 
to the y axis, therefore it only has sensitivity to the m g — m b difference, but little sensitivity to the absolute value of 
the gluino mass. The size of the band with significant probability is also different event by event, which means that 
some events will have more weight in the determination of the mass parameters. 

For each event, one defines in this way curves of equal probability in the (nig,m b ) plane. By combining the 
probabilities for different events, a region of maximum probability in the (nig,mt) is found, where the curves of 
maximum probability for all events approximately cross. Given the fact that the curves have different shapes for 
different events, the region thus defined has a limited size, for an adequate number of events. This region can be taken 
as a measurement for nig and m b . We perform the combination as the product of C, for the all events, and we define 

log£ C omb(mg,mj) = \og£(m g ,m~ b ) (13) 

events 

As an estimator of the probability for a given (m,g,m b ) pair, one can use 

^X 2 ( m ~gi m l) ~ Alog£ ee log£ max - \og£ comh (m g ,m b ) (14) 
Given the approximations introduced this does not however correspond to the statistical definition of A% 2 - 



B. Event analysis 

By following the procedure described in the previous subsection we can build the combined likelihood for all the 
events defined as: 

log C comh (m g , m g + Am g ,m b ,m b + Am b ) = ^ log(£(m§, m~ g + Am s , mg, m~ h + Amg) + c). (15) 

events 

As in Eq. (O we have introduced a constant cutoff parameter c = 0.001 and the number of Monte Carlo experiments 
used to build the likelihood for each event n — 10000. 

As we can see in Fig.|2 there are significant backgrounds from accidental leptons. The background subtraction must 
be carried out using events with OSOF lepton pairs. The correct log likelihood function is schematically expressed as 

log£ tota i(m§,mg, )= log£(m sig ,m bg ) + ^ £ ( m bg) (16) 

OSSF OSOF 
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where m s ; g express the parameters relevant to the signal distribution, such as the masses of the sparticles involved in 
the cascade decay, the decay branching ratios and so on. On the other hand mb g are all the other parameters relevant 
to the OSSF and OSOF events. This is rather complex procedure which is out of the scope of this paper. Instead, 
we take the difference of the functions for OSSF and OSOF lepton pair events 

log£ su b = log£ ssF - log£osoF = ^2 log£- ^2 log£. (17) 

OSSF OSOF 

In the limit of infinite statistics, log£ su b should be independent from the contribution of accidental lepton pairs. 
Therefore we use log £ S ub in this paper. 

We plot the contours of the function log £ su b in Fig. where plots (a) and (b) [(c) and (d)] are for tan /3 = 10 
[tan/3 = 20]. The distributions (a) and (c) are produced accepting all the events which pass the selections, whereas 
distributions (b) and (d) are produced using an event sample where the events including a 62 decay have been rejected. 

In Fig. [5] (a) and (c), the position of the peak for m g — m~ b is roughly consistent with the input value. Unlike 
the gluino and sbottom mass fits in the previous section, we obtain the correct peak position without the need of 
artificially choosing among multiple solutions. The likelihood distribution can be used to determine the g and b. We 
restrict the likelihood distribution for 591 GeV < m g < 599 GeV(within 4 GeV from the input gluino mass). We then 
fit the distribution around the peak assuming gaussian distribution, The likelihood distribution peaks at the gluino 
and sbottom mass difference as 99.5 GeV for tan/3 = 10, 104.2 GeV for tan/3 = 15, and 113.9 GeV for tan/3 = 20, 
where the input value is 103.3 GeV, 109.9 GeV and 116.5 GeV, respectively. The fitted values display shift of about 
4 GeV from the true value. We ascribe this effect to our simplified modeling of the jet smearing in building the 
likelihood function, which should disappear once the detector response is properly taken into account in the unfolding 
procedure. 

By comparing the left side with the right side of Fig [SI we also observe a slight shift in the position of the maxima 
of the distributions, showing that the distributions are sensitive to the presence of 62 decays. The 62 contribution 
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(a) tan/S=10, all (b) tan/S=10, no sbottom2 




Mass difference (GeV) Mass difference (GeV) 



(c) tan/S=20, all (d) tan/3=20, no sbottom2 




Mass difference (GeV) Mass difference (GeV) 

FIG. 5: Contours of the likelihood function log£ au b in the (m§ —mt,mg) plane: (a) and (b) for tan/3 = 10 and (c) and (d) for 
tan/3 = 20, respectively. The contours (b) and (d) are made without 62 contributions. 

however manifests itself in Fig. [5] (a) and (c) only as a flattening of the distribution around m§ — mz = 70 GeV at 
nig ~ 595 GeV. No secondary peak can be observed because of the experimental smearing, and of the fact that the 
branching ratio into 62 is much smaller. In Fig. E{a), we show the distribution of log£ su b as a function of rag — 
at tan/3 = 20, restricting the gluino mass in the region 591 GeV< m§ <599 GeV again. On the left of the peak 
corresponding to the b\ mass, we see a small bump in the distribution. This bump is not observed in the mass 
distribution made without 62 contribution (Fig. Etb)). In order to claim the presence of a second component in the 
distribution on the data, the ability of correctly reproducing the likelihood distribution for b\ events would be needed. 
It is also difficult to extract a statistical significance for the bi shoulder as our definition of the likelihood function is 
approximate one, and we did not treated the background subtraction correctly as can be seen in Eq. I|17|) . 

In Fig. |2| it is rather hard to see the effect of 62 unlike in Fig.El The apparent discrepancy probably comes from the 
fact that the likelihood analysis is more sensitive to the model parameter than simply solving Eq. Q for a fixed gluino 
mass. The likelihood analysis not only takes care of the most plausible value of the sbottom mass for a fixed gluino 
mass, but also includes possible statistical fluctuations, which vary event by event as seen in Fig. 21 For example, 
badly mis-measured events have less chance to be consistent with the input gluino mass, providing a natural cut for 
the event selection. 

A significant part of the background under the 62 is the tail of the smeared probability distribution for events which 
correctly reconstruct the b\ mass. One can therefore try to remove from the distributions the events consistent with 
b\ in order to improve the signal to background ratio for 62- To reduce the b\ events, an event is required to satisfy 
the condition 

C< £ ™t (18) 

cut region 

where the sum is made for bins in a cut region in the (m§ — mg,m§) plane. We choose the region as 
550 GeV< rrig < 650 GeV, and m(min) < m g — < m(max), which corresponds to the region around the 61 
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FIG. 6: The likelihood as a function of m§ — for tan/3 = 20: (a) for all events and (b) without 62 events. 





m(min) 


m(max) 


all 


62 only 


tan j3 = 


10 


87 


108 


64.1 ± 0.2 


62.7 ± 0.3 






82 


113 


60.3 ± 0.3 


60.2 ± 0.5 


tan (3 = 


15 


94 


115 


61.1 ± 0.7 


65.3 ± 0.4 






89 


120 


62.7 ± 0.3 


62.6 ± 0.6 



TABLE V: The peak positions at the smaller m,g — region in Fig. [S]for tan — 10 and 15. The peak positions for the signal 
62 distribution are also shown. 



peak. We use the cut value £ cu t = 20. The relevant m(min) and m(max) values are listed in Table Ivl The contours 
of log£ su b after this cut are shown in Fig.0 (a) and (b) for tan/3 = 10 and (c) and (d) for tan/3 = 15. The contours 
(b) and (d) are made without the 6 2 contribution. By comparing the contours with and without the b 2 contribution, 
the presence of the 62 can be clearly observed in the plots. In Fig. [SJ we again plot the likelihood function as the 
function of m g — mg, for 550 GeV< m g < 650 GeV. The distributions without the 6 2 contribution are also shown as 
dashed histograms. The signal to background ratio is much improved with respect to what is seen on the left side of 
Fig. El and it is about 1:1. However, from an inspection of the mass distribution event by event, the purity of the 
signal after the likelihood cut does not appear significantly improved with respect to Fig. UJ Moreover the position 
of the peak corresponding to 62 is dependent from the cuts applied. This is illustrated in Table [V] where the results 
of fits to the peak position are shown for two different values of the applied cuts both for the full sample and for the 
pure 62 signal. Even when only considering the 62 signal, the peak position depends on the cuts, albeit with a milder 
dependence than for the full sample. It will therefore be problematic to extract a mass measurement from Fig. [SJ 
even assuming a priori the existence of a 62 contribution. 



VII. DISCUSSION 



Supersymmetric models predict the existence of heavy superpartners which decay subsequently into the lighter 
superpartners. The lightest SUSY particle is stable and neutral, and escapes detection. Therefore two undetected 
particles will be present in each event. Moreover the partonic center of mass energy is unknown for hadron collisions. 
As a result, the complete kinematic reconstruction of SUSY events at hadron colliders is problematic. 

We propose a new analysis method to solve the decay kinematics of this decay at the LHC, based on imposing the 
on-shell condition on the momenta of the particles participating in the cascade. The single cascade decay is solvable if 
a decay chain consisting of at least 4 successive two-body decays, involving 5 sparticles can be identified. In this case 
each event defines the 4 dim hypersurface in the 5 dim sparticle mass space. The potential of this method should be 
compared to the method previously used for this analysis based on the measurement of kinematic edges of invariant 
mass combinations of the detected decay particles. 

The merits and demerits of the new method may be summarized as follows: 
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(a) tan/S=10, all (b) tan/S=10, no sbottom2 




Mass difference (GeV) Mass difference (GeV) 



(c) tan/S=15, all (d) tan/S=15, no sbottom2 




Mass difference (GeV) Mass difference (GeV) 

FIG. 7: The contours of the likelihood function in the (nig — mg, m§) plane after the cut given in Eq. 1)18^ : (a) and (b) are 
with and without 62 events at tan/3 = 10, respectively, and (c) and (d) at tan/3 = 15. 

• The cascade decay is solved based on the exact formulae. 

• Mass peaks would be reconstructed, as opposed to the kinematic edges. This allows us to perform measurements 
even if the significant backgrounds exist or statistics are small. Note that backgrounds will not exhibit peaks 
corresponding to the signal region. 

• In the case where sfermion masses are heavier than gaugino masses, the cascade decays are expected to be 
shorter, and one can not therefore use this method. The endpoint method provides a mass information even in 
this case. 

We note that in case two decay chains can be simultaneously identified in the event, the method can be applied 
to shorter decay chains, consisting of only three decays each. In fact in this case two additional constraints can be 
applied by requiring that the sums of the transverse components of the momenta of the two lightest neutralinos equal 
the two components of the measured missing transverse momentum. 

Further on, if one sparticle cascade decay is solved by the mass relation method, the candidate LSP momentum px 
would be obtained. We can then calculate the transverse momentum of the other LSP p' T as: 

FT = ~PT + -Pmiss- (19) 

For the cascade decay to which the second LSP belongs only two components of the neutralino four momentum are 
unknown, therefore a cascade decay with n^ ecay > 3 can be solved. 

To see the performance of the mass relation method, we have studied in this paper the problem of measuring 
sparticle masses in the cascade decay: g — > 662 — * — ► lb\b%l% — > Xi^i^2-^i^2- We have performed a detailed 
study for some benchmark SUSY model points including backgrounds and a parameterized simulation of detector 
effects. We performed the exercise in a simplified fashion, by fixing the masses of the three lighter particles to avoid 
the practical complications in handling a large number of parameters. Then each event becomes an allowed curve in 
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tan/?=10 



tan/?=15 




60 80 100 120 140 

Mass difference (GeV) 




60 80 100 120 140 

Mass difference (GeV) 



FIG. 



■mi with the cut given in Eq. I|18jl : (left) for tan/3 = 10 and (right) for tan/3 = 15. 



the gluino-sbottom mass plane, passing through the point corresponding to the true gluino and sbottom masses. We 
first addressed the mass reconstruction through the "event pair analysis" , which determines the sparticle masses from 
the distribution of the solutions of any event pairs in the selected sample. The method reconstructs correctly 
for SPSla where tan/3 is varied from 10 to 20. Note that the signal branching ratio becomes a factor 4 smaller for 
tan/3 = 20 with respect to tan/3 = 10, but the b\ peak is still clearly observable. On the other hand, in order to 
obtain the correct mass, one needs to artificially choose among the multiple available solutions. 

A more global approach requires the usage of all available events simultaneously. For this approach we constructed 
an approximate likelihood function for the true gluino and sbottom masses, taking into account the experimental 
smearing in the measurement of the 6-jets. When the peak position of the likelihood distribution is used to extract 
the mass, the b\ mass is measured without the problem of the multiple solutions of the event pair analysis. 

We also try to probe the presence of a 62 in our Monte Carlo sample. The sbottom mass matrix is parameterized 
by three parameters mg , and the mixing angle 9. Successful extraction of the 62 would be an important step 

to fully understand the nature of the third generation sparticles. No clear result is achieved for 62, as the branching 
ratio into 62 is small for the parameters we have chosen, and also the small difference between the two sbottom states 
is comparable to the resolution in the experimental measurement of jet momenta. It is however clear that, even with 
the small statistics available for the case tan/3 = 20, a hint for the deviation from a single-mass case can be seen in 
the distribution. 

The physics output of our analysis is therefore the possibility to extract information on the third generation 
sector, even for rather small input statistics. The measurement of the third generation sparticle masses are important 
theoretically. In mSUGRA, &l(_r) masses are same as the other sparticle masses at the GUT scale but non-universality 
is induced by RGE running at the weak scale due to the Yukawa coupling. In addition to that, the left-right mixing 
of sbottom is induced by the F term of the Superpotential which is proportional to /xtan/3. The mass shift around 
20 GeV from tan/3 = 10 and tan/3 = 20 is due to the mixing effect. Our method is sensitive to the tan/3 dependence 
as it could yield a sensitivity to the mass difference between the two states of the order of a few GeV. 
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As described in section 2, we have solved the cascade 
lomentum bv the momenta of U . lo and Jh as in Eq. J3J 



Appendix 

decay in Eq. JIJ by expanding the the lightest neutralino 



as described m section z, we nave solved tnc 
momentum by the momenta of Zj, Z 2 and &i as in 

P; 

The parameters (a, 6, c) may be written as the function of sparticle masses and the LSP energy as 



p^o = ap ei + bpt 2 + cp bl 




where 



[xi]m|o + [x 2 ]m? + [x 3 ]m x o^o 



X = [xi,x 2 ,x 3 ], 

/ — m~ + to? 

1 



= M~ l X 



2 
b 



X.1 Ai 



[xi] = 



ix 2 j = 



and M is defined already in section 2 as: 

M 




m% + to? + 2p £l ■ p e 



a.2 <■ 

[x 3 ] = £ 



TO x? 



Ki ' Ki Wi ' W 2 Wi • Pfci 

= I Kr P£ 2 Pe 2 ■ Pi 2 Pe 2 ■ 



, P^i ' Pbi Pt 2 ■ Pfci Pbi • Pbi 

By using the on-shell condition of the neutralino mass 



E% = (a, 6, c)M 




we 



obtain the following equation: 



A, 



^x? 



■21 A 



13 "I T - A23 

XY / 



^x? 



Xi 



^ TO? 



m 



to? 
x? 



4 
xY 



6 ^22 



where A.y = [x;] T .A/f 1 [xj] — 6ij(i — 2)/m? The solution of E%a is expressed 



as 



^x? = 



TO-r^O 
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D = [A 



13 
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Xi 
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TO x? 



rn 



x? 



f ^22 



(4| 3 _ A 33 A 22 ) + 2 {A 13 A 23 - A 33 A 12 ) 

V m -^o / 



0, 



When the sbottom comes from gluino decay, we can further use the gluino mass shell condition: 

m j = 2 Pb 2 ■ Px° + 2 Pb 2 ■ (Pbi +Pt 2 +Pil) + m \ + m h 



(20) 



(21) 



(22) 



(23) 
(24) 
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where 
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The equation involving the gluino and sbottom masses is of the form 

Qiim| + 2<9i2m|mg + Q 22 m^ + 2Qim~ + 2Q 2 to? + Q = 0, 

where 



111 



1, Q 12 = -2F 1 -1, 



(25) 



n - rtP i 1\2 p 2 ^23 ~ ^33^22 

<y 22 - {2*1 + 1) - * 2 j2 , 

^33 

Qi = -2F to|o - 2p b2 ■ (p bl + p tl +pe 2 ) - ml, 

r\ r\ r\ r 2 2 ^13^23 ~ ^33^12 

1 A33 

Q = ^-^4 ^"i lA - (28) 

1 ^33 

It should be noted that one can derive a equation of the form m| + C\m\ + C2 — for any particle / involved in the 
cascade. 
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